<html>
<head>
  <meta HTTP-EQUIV="Content-Type" CONTENT="text/html;charset=ISO-8859-1">
  <title>rbfpreimg2.m</title>
<link rel="stylesheet" type="text/css" href="../../../m-syntax.css">
</head>
<body>
<code>
<span class=defun_kw>function</span>&nbsp;<span class=defun_out>z&nbsp;</span>=&nbsp;<span class=defun_name>rbfpreimg2</span>(<span class=defun_in>varargin</span>)<br>
<span class=h1>%&nbsp;RBFPREIMG2&nbsp;RBF&nbsp;pre-image&nbsp;problem&nbsp;by&nbsp;Gradient&nbsp;optimization.</span><br>
<span class=help>%</span><br>
<span class=help>%&nbsp;<span class=help_field>Synopsis:</span></span><br>
<span class=help>%&nbsp;&nbsp;z&nbsp;=&nbsp;rbfpreimg2(model)</span><br>
<span class=help>%&nbsp;&nbsp;z&nbsp;=&nbsp;rbfpreimg2(model,options)</span><br>
<span class=help>%</span><br>
<span class=help>%&nbsp;&nbsp;<span class=help_field>Description:</span></span><br>
<span class=help>%&nbsp;&nbsp;&nbsp;z&nbsp;=&nbsp;rbfpreimg2(model)&nbsp;it&nbsp;uses&nbsp;gradient&nbsp;method&nbsp;to&nbsp;solve&nbsp;</span><br>
<span class=help>%&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;the&nbsp;pre-image&nbsp;problem&nbsp;for&nbsp;the&nbsp;Radial&nbsp;Basis&nbsp;Function&nbsp;(RBF)&nbsp;</span><br>
<span class=help>%&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;kernel.&nbsp;The&nbsp;function&nbsp;'fminunc'&nbsp;of&nbsp;the&nbsp;Matlab&nbsp;Optimization&nbsp;</span><br>
<span class=help>%&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;toolbox&nbsp;is&nbsp;exploited&nbsp;for&nbsp;1D&nbsp;search&nbsp;along&nbsp;the&nbsp;gradient&nbsp;</span><br>
<span class=help>%&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;direction.</span><br>
<span class=help>%</span><br>
<span class=help>%&nbsp;&nbsp;&nbsp;z&nbsp;=&nbsp;rbfpreimg2(model,options)&nbsp;use&nbsp;to&nbsp;specify&nbsp;the&nbsp;control</span><br>
<span class=help>%&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;parameters&nbsp;of&nbsp;the&nbsp;gradient&nbsp;optimization.</span><br>
<span class=help>%&nbsp;</span><br>
<span class=help>%&nbsp;<span class=help_field>Input:</span></span><br>
<span class=help>%&nbsp;&nbsp;model&nbsp;[struct]&nbsp;Kernel&nbsp;expansion:</span><br>
<span class=help>%&nbsp;&nbsp;&nbsp;.Alpha&nbsp;[num_data&nbsp;x&nbsp;1]&nbsp;Weight&nbsp;vector.</span><br>
<span class=help>%&nbsp;&nbsp;&nbsp;.sv.X&nbsp;[dim&nbsp;x&nbsp;num_data]&nbsp;Vectors&nbsp;determining&nbsp;the&nbsp;kernel&nbsp;expansion.</span><br>
<span class=help>%&nbsp;&nbsp;&nbsp;.options.arg&nbsp;[1x1]&nbsp;Argument&nbsp;of&nbsp;the&nbsp;RBF&nbsp;kernel&nbsp;(see&nbsp;'help&nbsp;kernel').</span><br>
<span class=help>%</span><br>
<span class=help>%&nbsp;&nbsp;options&nbsp;[struct]&nbsp;Control&nbsp;parameters:</span><br>
<span class=help>%&nbsp;&nbsp;&nbsp;.min_improvement&nbsp;[1x1]&nbsp;Minimal&nbsp;allowed&nbsp;improvement&nbsp;of&nbsp;the&nbsp;objective&nbsp;</span><br>
<span class=help>%&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;function&nbsp;in&nbsp;two&nbsp;consecutive&nbsp;steps&nbsp;(default&nbsp;1e-3).</span><br>
<span class=help>%&nbsp;&nbsp;options.start_t&nbsp;[1x1]&nbsp;Starting&nbsp;value&nbsp;of&nbsp;the&nbsp;1D&nbsp;search&nbsp;procedure&nbsp;</span><br>
<span class=help>%&nbsp;&nbsp;&nbsp;(default&nbsp;1e-3).</span><br>
<span class=help>%</span><br>
<span class=help>%&nbsp;<span class=help_field>Output:</span></span><br>
<span class=help>%&nbsp;&nbsp;z&nbsp;[dim&nbsp;x&nbsp;1]&nbsp;Found&nbsp;preimage.</span><br>
<span class=help>%</span><br>
<span class=help>%&nbsp;See&nbsp;also&nbsp;</span><br>
<span class=help>%&nbsp;&nbsp;RBFPREIMG,&nbsp;RBFPREIMG3,&nbsp;RSRBF,&nbsp;KPCAREC.</span><br>
<span class=help>%</span><br>
<hr>
<span class=help1>%&nbsp;<span class=help1_field>About:</span>&nbsp;Statistical&nbsp;Pattern&nbsp;Recognition&nbsp;Toolbox</span><br>
<span class=help1>%&nbsp;(C)&nbsp;1999-2003,&nbsp;Written&nbsp;by&nbsp;Vojtech&nbsp;Franc&nbsp;and&nbsp;Vaclav&nbsp;Hlavac</span><br>
<span class=help1>%&nbsp;&lt;a&nbsp;href="http://www.cvut.cz"&gt;Czech&nbsp;Technical&nbsp;University&nbsp;Prague&lt;/a&gt;</span><br>
<span class=help1>%&nbsp;&lt;a&nbsp;href="http://www.feld.cvut.cz"&gt;Faculty&nbsp;of&nbsp;Electrical&nbsp;Engineering&lt;/a&gt;</span><br>
<span class=help1>%&nbsp;&lt;a&nbsp;href="http://cmp.felk.cvut.cz"&gt;Center&nbsp;for&nbsp;Machine&nbsp;Perception&lt;/a&gt;</span><br>
<br>
<span class=help1>%&nbsp;<span class=help1_field>Modifications:</span></span><br>
<span class=help1>%&nbsp;10-jun-2004,&nbsp;VF</span><br>
<span class=help1>%&nbsp;03-dec-2003,&nbsp;VF</span><br>
<br>
<hr>
<span class=comment>%&nbsp;process&nbsp;input&nbsp;arguments</span><br>
<span class=comment>%----------------------------------</span><br>
<span class=keyword>if</span>&nbsp;<span class=stack>nargin</span>&nbsp;&gt;&nbsp;2,<br>
&nbsp;&nbsp;z=foo(<span class=stack>varargin</span>{1},<span class=stack>varargin</span>{2},<span class=stack>varargin</span>{3},<span class=stack>varargin</span>{4},<span class=stack>varargin</span>{5});<br>
&nbsp;&nbsp;<span class=jump>return</span>;<br>
<span class=keyword>else</span><br>
&nbsp;&nbsp;model&nbsp;=&nbsp;<span class=stack>varargin</span>{1};<br>
&nbsp;&nbsp;<span class=keyword>if</span>&nbsp;<span class=stack>nargin</span>&nbsp;&lt;&nbsp;2,&nbsp;options=[];&nbsp;<span class=keyword>else</span>&nbsp;options&nbsp;=&nbsp;c2s(&nbsp;<span class=stack>varargin</span>{2});&nbsp;<span class=keyword>end</span><br>
&nbsp;&nbsp;<span class=keyword>if</span>&nbsp;~isfield(options,<span class=quotes>'min_improvement'</span>),&nbsp;options.min_improvement&nbsp;=&nbsp;1e-3;&nbsp;<span class=keyword>end</span><br>
&nbsp;&nbsp;<span class=keyword>if</span>&nbsp;~isfield(options,<span class=quotes>'start_t'</span>),&nbsp;options.start_t&nbsp;=&nbsp;1e-3;&nbsp;<span class=keyword>end</span><br>
&nbsp;&nbsp;<span class=keyword>if</span>&nbsp;~isfield(options,<span class=quotes>'attempts'</span>),&nbsp;options.attempts&nbsp;=&nbsp;10;&nbsp;<span class=keyword>end</span><br>
<span class=keyword>end</span><br>
<br>
[dim,num_sv]=size(model.sv.X);<br>
ker&nbsp;=&nbsp;<span class=quotes>'rbf'</span>;<br>
arg&nbsp;=&nbsp;model.options.arg;<br>
iXi&nbsp;=&nbsp;sum(&nbsp;model.sv.X.^2)';<br>
s2&nbsp;=&nbsp;arg^2;<br>
<br>
<span class=comment>%&nbsp;Selection&nbsp;of&nbsp;the&nbsp;starting&nbsp;point&nbsp;out&nbsp;of&nbsp;the&nbsp;model.sv.X.</span><br>
<span class=comment>%&nbsp;The&nbsp;point&nbsp;in&nbsp;which&nbsp;is&nbsp;the&nbsp;objective&nbsp;function&nbsp;minimal&nbsp;is&nbsp;taken.</span><br>
<span class=comment>%&nbsp;Minimum&nbsp;over&nbsp;50&nbsp;randomly&nbsp;drawn&nbsp;points&nbsp;is&nbsp;used.</span><br>
<span class=comment>%--------------------------------------------------------------</span><br>
<br>
rand_inx&nbsp;=&nbsp;randperm(&nbsp;num_sv&nbsp;);<br>
rand_inx&nbsp;=&nbsp;rand_inx(1:min([num_sv,50]));<br>
Z&nbsp;=&nbsp;model.sv.X(:,rand_inx);<br>
<br>
fval&nbsp;=&nbsp;kernel(Z,model.sv.X,ker,arg)*model.Alpha(:);<br>
fval&nbsp;=&nbsp;-fval.^2;<br>
<br>
[dummy,&nbsp;inx&nbsp;]&nbsp;=&nbsp;min(&nbsp;fval&nbsp;);<br>
z&nbsp;=&nbsp;Z(:,inx&nbsp;);<br>
<br>
<span class=comment>%&nbsp;Gradient&nbsp;descent&nbsp;optimization&nbsp;&nbsp;</span><br>
<span class=comment>%--------------------------------------</span><br>
change=inf;<br>
opt=optimset(<span class=quotes>'display'</span>,<span class=quotes>'off'</span>,<span class=quotes>'Diagnostics'</span>,<span class=quotes>'off'</span>,<span class=quotes>'LargeScale'</span>,<span class=quotes>'off'</span>);<br>
<span class=error>warning</span>&nbsp;off&nbsp;MATLAB:divideByZero;<br>
<span class=keyword>while</span>&nbsp;change&nbsp;&gt;&nbsp;options.min_improvement,<br>
&nbsp;&nbsp;&nbsp;<br>
&nbsp;&nbsp;&nbsp;<span class=comment>%&nbsp;compute&nbsp;gradient</span><br>
&nbsp;&nbsp;&nbsp;dotp&nbsp;=&nbsp;kernel(&nbsp;model.sv.X,z,ker,arg&nbsp;).*model.Alpha(:);<br>
&nbsp;&nbsp;&nbsp;dz&nbsp;=&nbsp;z*sum(dotp)&nbsp;-&nbsp;model.sv.X*dotp;<br>
<br>
&nbsp;&nbsp;&nbsp;<span class=comment>%&nbsp;auxiciliary&nbsp;variables</span><br>
&nbsp;&nbsp;&nbsp;zXi&nbsp;=&nbsp;model.sv.X'*z;<br>
&nbsp;&nbsp;&nbsp;dzXi&nbsp;=&nbsp;model.sv.X'*dz;<br>
&nbsp;&nbsp;&nbsp;<br>
&nbsp;&nbsp;&nbsp;Ai&nbsp;=&nbsp;-(1/(2*s2))&nbsp;*&nbsp;(iXi&nbsp;-&nbsp;2*zXi&nbsp;+&nbsp;z'*z);<br>
&nbsp;&nbsp;&nbsp;Bi&nbsp;=&nbsp;-(1/s2)&nbsp;*&nbsp;(z'*dz&nbsp;-&nbsp;dzXi);<br>
&nbsp;&nbsp;&nbsp;C&nbsp;=&nbsp;-(1/(2*s2))&nbsp;*&nbsp;(dz'*dz);<br>
&nbsp;&nbsp;&nbsp;<br>
&nbsp;&nbsp;&nbsp;<span class=comment>%&nbsp;1D-search&nbsp;to&nbsp;determine&nbsp;the&nbsp;size&nbsp;of&nbsp;the&nbsp;step&nbsp;in&nbsp;the&nbsp;gradient&nbsp;direction</span><br>
&nbsp;&nbsp;&nbsp;[t,fval]&nbsp;=&nbsp;fminunc(<span class=quotes>'rbfpreimg2'</span>,options.start_t,opt,model.Alpha,Ai,Bi,C);<br>
&nbsp;&nbsp;&nbsp;<br>
&nbsp;&nbsp;&nbsp;old_z&nbsp;=&nbsp;z;<br>
&nbsp;&nbsp;&nbsp;z&nbsp;=&nbsp;z&nbsp;+&nbsp;dz*t;<br>
&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;<br>
&nbsp;&nbsp;&nbsp;change&nbsp;=&nbsp;sum((z-old_z).^2);<br>
<span class=keyword>end</span><br>
<br>
<span class=jump>return</span>;<br>
<br>
<span class=comment>%---------------------------------</span><br>
<span class=defun_kw>function</span>&nbsp;<span class=defun_out>f</span>=<span class=defun_name>foo</span>(<span class=defun_in>t,Alpha,Ai,Bi,C</span>)<br>
<span class=comment>%&nbsp;</span><br>
<br>
f&nbsp;=&nbsp;Alpha(:)'*exp(Ai&nbsp;+&nbsp;Bi*t&nbsp;+&nbsp;C*t^2);<br>
f&nbsp;=&nbsp;-f^2;<br>
<br>
<span class=jump>return</span>;<br>
<span class=comment>%&nbsp;EOF</span><br>
<br>
<br>
</code>
